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ABSTRACT 

We have performed an analysis of the diffuse gamma-ray emission with the Fermi Large Area Telescope 
in the Milky Way Halo region searching for a signal from dark matter annihilation or decay. In the absence 
of a robust dark matter signal, constraints are presented. We consider both gamma rays produced directly in 
the dark matter annihilation/decay and produced by inverse Compton scattering of the e^/e~ produced in the 
annihilatioii/decay. Conservative limits are derived requiring that the dark matter signal does not exceed the 
observed diffuse gamma-ray emission. A second set of more stringent limits is derived based on modeling the 
foreground astrophysical diffuse emission using the GALPROP code. Uncertainties in the height of the diffusive 
cosmic -ray halo, the distribution of the cosmic -ray sources in the Galaxy, the index of the injection cosmic- 
ray electron spectrum and the column density of the interstellar gas are taken into account using a profile 
likelihood formalism, while the parameters governing the cosmic-ray propagation have been derived from fits 
to local cosmic -ray data. The resulting limits impact the range of particle masses over which dark matter 
thermal production in the early Universe is possible, and challenge the interpretation of the PAMELA/Fermj- 
LAT cosmic ray anomalies as annihilation of dark matter. 
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1. INTRODUCTION 

The nature of Dark Matter (DM) and its properties are 
still unknown, despite being one of the most widely inves- 
tigated topics in contemporary fundamental physics. How- 
ever, current and near future experiments are probing more 
and more of the parameter space predicted for the most pop- 
ular type of DM candidates. Weakly Interacting Massive Par- 
ticles ( WIMPs) (for a re view see Bergstrom| |2000[ |Bertone| 
|et al.|2004 ; Bertone|2010| ). In particular, high-energy gamma- 
ray astronomy can be used to search for signatures of DM 
in the Milky Way and beyond. Gamma-rays are products of 
hadronization and radiative loss processes, and are therefore 
unavoidably emitted in annihilation and decay of WIMPs. 
The propagation of gamma rays is mostly unaffected by the 
interstellar medium and Galactic magnetic fields, and there- 
fore the data retain information on the morphology of the 
emission region (unlike, e.g., charged cosmic rays). The 
Large Area Telescope (LAT), onboard the Fermi gamma-ray 
observatory (Atwood et al. 2009 ), is now providing unprece- 
dented high quality gamma-ray data. 

We focus here on DM signatures in the diffuse gamma- 
ray emission as measured by Ferm/-LAT. About 90% of the 
LAT photons are of diffuse origin. The Galactic component 
encodes information on the propagation and origin of cos- 
mic rays, distribution of cosmic -ray sources, the interstellar 
medium, magnetic and radiation fields in our Galaxy and un- 
resolved point sources, while its extragalactic component pro- 
vides a signature of energetic phenomena on cosmological 
scales. Both components are expected to include a contribu- 
tion from DM annihilation/decay: the Galactic signal arises 
from the smooth DM halo around the Galactic Center and 
Galactic substructures, while the extragalactic one arises from 
the signal of DM annihilation processes throughout the Uni- 
verse integrated over all reds hifts. The extragal actic com- 
ponent is analyzed elsewhere ( Abdo et al.||2010a] ). Here we 
focus on searching for a potential signal from DM annihila- 
tion/decay in the halo o f the Milky W ay (for previou s work 
related to this topic see [Zhan g et al. 1(2009; Papuc ci & Stru- 
mia 2010; CirelH et al."2010; Malyshev et al. 201 If lBaxter"^ 
Dodelson.2011, .Lin et al. 2010, .Abbasi et al.,20llll. 



Due to the bright sources present in the Galactic Center and 
the bright di ffuse emissio n along the plane, it has been ar- 
gued, e.g. in Serpico & Zaharijas (2008), that the region of 



the inner Galaxy, extending 10°-20° away from the Galactic 
Plane, is promising in terms of the signal-to-background ratio 
S/N. As an additional advantage, the constraints on the DM 
signal in that region become less sensitive to the unknown 
profile of the DM halo. In particular the S/N of cored pro- 
files is only a factor ~ 2 weaker tha n for the Navarro-Frenk- 
White (NFW) profile ( [Navarro et al.,199 6 ) in this region. This 
should be compared to an order of magnitude of uncertainty 
when one considers the Galactic Center region. The increase 
in S/N ratio away from the plane is further emphasized for 
DM models in which DM annihilations result in a significant 
fraction of leptons in the final state. These leptons propagate 
in the Galaxy and produce high-energy gamma rays mainly 
through inverse Compton scattering on the interstellar radia- 
tion field. By diffusing away from the Galactic Center region, 
electrons produce an extended gamma-ray signal which fur- 
ther enhan ces the S/N at higher Galactic latitudes ( Borriello 
[et al.|2 009b). 

Consequently, we investigate here a large region of interest 
covering the central part of the Galactic Halo, while masking 



out the Galactic Plane. We test the LAT data for a contribution 
from the DM annihilation/decay signal by a fit of the spectral 
and spatial distributions of the expected photons in the region 
of interest to the LAT data. In doing so, we take into account 
the most up-to-date mod eling of the diff use signal of astro- 
physical origin (Ackerm ann et al.[[2012a ), adapting it to the 
problem in question. This paper is organized in the follow- 
ing way: In Sec. j2]we describe our modeling of the diffuse 
gamma-ray emission and the way we parameterize it. Sec. |3] 
outlines our general approach to fitting for DM signals in the 
presence of uncertainties in the astrophysical foregrounds. In 
Sec. |4]we describe the DM and gas maps used in this work 
while in Sec. |5|we define our dataset and region of interest. 
In Secl6]we derive DM limits without modeling of the back- 
grouncrwhile Sec. |7] contains the details of the fit procedure 
for the limits which include modeling of the background. The 
results are presented in Sec. [8] while Sec. [9] contains further 
discussions on the background model uncertainties. In Sec. 
[TOlwe summarize and conclude. 

2. MODELING OF THE HIGH-ENERGY GALACTIC 
DIFFUSE EMISSION 

The Galactic diffuse gamma-ray emission is produced by 
the interaction of energetic cosmic-ray (CR) electrons and nu- 
cleons with the interstellar gas and radiation field. Its main 
components are photons from the decay of neutral pions pro- 
duced in the interaction of the CR nucleons with the interstel- 
lar gas, bremsstrahlung of the CR electron (CRE) population 
on the interstellar gas and their inverse Compton (IC) scatter- 
ing off the Interstellar Radiation Field (ISRF). Efficient mod- 
eling of the diffuse gamma-ray emission needs an accurate 
description of both the interstellar gas and radiation targets as 
well as the distribution of CRs in the Galaxy. 

The gas in the Interstellar Medium (ISM) consists mostly 
of atomic hydrogen (H l) and, to a lesser extent, molecu- 
lar hydrogen (H2) which is concentrated along the Galac- 
tic Plane. Ionized (H ll) hydroge n is subdominant , although 
it has a larger vertical sc ale (see Ackermann et al. (2012a); 
[Moskalenko et al.] ( [2004| l for more details). Helium (He) is 
also important, being ~ 25% by mass of the ISM, and its 
distribution is assumed to follow that of interstellar hydrogen 
( [Ackermann et al.|20I2a| l. Velocity resolved radio surveys of 
the 21 cmhyperfine structure transition of H I and correspond- 
ing surveys of the 2.6 mm CO J(l — > 0) transition (using the 
CO density as a proxy for the H2 density) are used to build 
maps of the interstellar gas in different annuli ( Ackermann] 



etal. 



2012a I providing effectively a 3 



J^model of the gas dis- 
tribution in the Galaxy. The conversion factors Xco between 
CO line intensity and H2 column density have been observed 
to vary throughout the Galaxy (Abdo et al. 2010d). Total gas 
colum n density es timated from E(B- V) visual reddening maps 
( |Schleg el et al.ll998 ) has been shown to be complement ary to 
the one estimated from H I and C O surveys combined ([Gre-| 
[nier et al.,2005f As described in Ackermann et al.[p012a[ ), 
we take this into account by correcting the gas column den- 
sity for each line of sight according to the value derived from 
the E(B- V) map, except in the regions of high extinction (see 
Sec.|i2|. 

We use a 2Dh-1 cylindrically symmetric model (2 spatial 
dimensions and the frequency dimension) of the ISRF, com- 
puted based on a model of the radiation emission of stel- 

More precisely the mode! is only pseudo 3D due to the near-far ambi- 
guity in the inner Galaxy I Ackermann et al. 2012a I. 
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Figure 1. Upper panel: Spatial (left) and spectral (right) distribution of gamma rays originating from tlie annihilation of a 250 GeV WIMP into bb. The 
left figure shows the expected intensity at E=10 GeV for the full sky in Galactic coordinates. A NFW profile is assumed for the DM halo and a value of 
{(T^v) = 4 X 10^-^cm^s"' for the DM annihilation cross section. For comparison purposes typical spectra of the astrophysical emission from tt" decay and 
Inverse Compton (IC) scattering are displayed in the right figure. The map also shows the boundaries of the region used to plot the average spectra of the right 
panel, and which we will use for the analysis described in this work. Central panel: Same for a 250 GeV WIMP annihilating into ^jT ■ The contribution from 
IC and from Final State Radiation (FSR) are shown separately in the spectrum and are superimposed in the spatial distribution. Lower panel: Spatial (left) and 
spectral (right) distribution of the IC emission of an astrophysical CR source population distributed uniformly in Galactocentric radius within 1 kpc from the 
Galactic Center and with a scale height of 200 pc. 



lar populat ions and furth er reprocessing in the Galactic dust 
([Moskalen ko et al.|2006] l. 

We use the GALPROP code ( [Strong et al |2000| l v54, to cal- 
culate the propagation and distribution of CRs in the Galaxy. 

The code is further used to create sky maps of the expected 
gamma-ray emission from the interactions of the CRs with 
the ISM and ISRF based on the models of the gas and radia- 
tion targets described above. GALPROP approximates the CR 
propagation by a diffusion process into a cylindrical diffusion 



zone of half-height z/, and radius /?/,. CREs and nuclei are in- 
jected by a parametrized distribution of CR sources. Energy 
losses, production of secondary particles in interactions and 
reaccelerati on of CRs in the I SM are taken into account (for 
details see Strong et al. '2000V Several important parameters 
enter the GALPROP modeling: the distribution of CR sources, 
the half-height of the diffusive halo Zh, the radial extent of 
the halo Rh, the nucleon and electron injection spectrum, the 
normalization of the diffusion coefficient Dq, the rigidity de- 
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pendence of the diffusion coefficient S (D{p) = (SDoip / pof^ 
with pq = 4 GV being the reference rigidity and j3 = v/c) the 
Alfven speed va (parametrizing the strength of re-acceleration 
of CRs in the ISM via Alfven waves) and the velocity of the 
Galacti c winds perpendicular to t he Galactic Plane Vc- Fol- 
lowing Ackermann et al. (20 12a I we will parameterize the 
nucleon mjection spectrum as a broken power law in rigidity 
with 7p,i,7p,2 the index of the spectrum before and after the 
break respectively and p^,, p the break rigidity. Similarly, the 
electron injection spectrum is parametrized as a double bro- 
ken power law with 7f, i , 7^ ,2 , 7f..3 , the index of the spectrum in 
the three rigidity zones and Pbr,e.\ , Pbr.,e,2 the two breaks. 



We use in the following the results and formalism of Ack- 



ermann et al. ( 2012a| l and we briefly summarize here the ap 



preach and resul ts presented there. The reader is referred to 
I Ackermann et a l. (2012a) for more details and thorough dis- 
cussion. In that work a grid of models is considered with 4 
values of Zh, 2 values of R/j and 4 different models of CR 
Source Distributions (CRSDs). The CRSDs are set to corre- 
spond to the incompletely determined distributions of Super- 
nova Remnants (SNR), o r tracers of star formation a n d col- 
lapse (pulsai's, OB stai's) Case & Bhattacharya| ( |1998| l; |Yusi-| 
[?ov & Ktigtik (2004). For each of these 32 models 4 different 
assumptions orTtHe column density of the hydrogen gas de- 
rived from its tracers are made for a total of 128 models. For 
each of the models, a fit of the model prediction to the local 
intensity of different CR nuclei and the B/C ratio is performed 
in order to fix the parameters (Do , va , jp. i , 7/,.2 , Phr)- Thus, the 
CR fit provides the injection spectrum for nuclei, the Alfven 
speed and the relation (z/,,Do) for different values of zi,- In a 
second step they determine also the electron injection spec- 
trum from a fit of the model to the local spectrum of CRE 
using the diffusion parameters obtained in the first fit. How- 
ever, we are not going to use the result of the second step, 
since we will instead fit the electron spectrum from gamma- 
ray data. As the last step in Ackermann et al. (2012a) an all- 
sky fit to the Fenni-LAI gamma-ray data is then performed 
to find the remaining parameters like the Xco factors^ The 
ISRF normalizations in various regions of the sky are also left 
as parameters free to vary in the fit, to account for possible 
uncertainties in the ISRF itself and the CRE distribution. In 
our analysis, we will also consider the ISRF uncertainties in 
more detail in Sec. |9] 

Thu s, we use only the results f rom the first step of the anal- 
ysis inl Ackermann et al.| ( 2012a| l described above, but then al- 
low for more freedom in certam parameters governing the CR 
distribution and astrophysical diffuse emission and constrain 
these parameters by fitting the m odels to the LA T gamma- 
ray data. Compared to Ackermann et al. (2012a|, the main 
difference in our analysis is the use of a free CRSD whose 
parameters will be determined from the fit t o the LAT data, as 
opposed to the 4 fixed CRSDs explored in [Ackermann et al.| 
( [20T2a ). The procedure which we employ to obtain the best 
lit CRSD is described in detail in Sec. 17] Another impor- 
tant difference is that we will keep the CRE source distribu- 
tion and proton source distribution (which we will refer to as 
eCRSD and pCRSD) separate. This is justified as we don't 
know a priori if the bulk of CR protons and electrons is in- 
jected by the same class of sources. We report in Table [T]the 
CR diffusion and injection parameters f or different values of 
Zh taken from [Ackermann et al.| ( |20 125) 1 which we will use in 

The procedure is iterated few times in order to consider the feedback of 
the renormalized Xco factors in the propagation of CREs. 



the foUowinj^ Note that the diffusion parameters in prin- 
ciple depend on the CRSD, but the dependence is weak and 
will be neglected in the following. There is also slight depen- 
dence on the parameters used to produce the gas maps (see 
Sec. [4| and the assumed Xco distribution, which is also weak 
and will be neglected as well. Besides the free CRSDs and 
the scan over different values of zi, we will also scan electron 
injection spectra by varying the index 7^ 2 wh ile we will fix 
7(,,i = 1.6, 7e,3 = 4 (Ackermann et al. 2012a| l, pbr,e.i = 2500 
MV and pbr,e,3 = 2.2 I V. ihe last two parameters are left free 



to vary in the analysis performed in [Ackermann et aL] ( 2012a| l 
although the fitted values typically differ by less than 20% for 
Pbr.e.i and less than 10% phr.e,3 with respect to the values we 
report above. 

As a sanity check it should be also verified a posteriori 
that the flux and spectrum of local protons obtained after the 
gamma-ray fit are consistent with the experimentally observed 
ones, since the pCRSD fit could, potentially, affect them, and 
this would make inconsistent the use of the initial diffusion 
parameters, based on the observed proton spectrum. We, in- 
deed, found that the proton spectrum after the second step is 
fully consistent with the observed one. We also verified that 
the normalization of the electron flux given by the gamma- 
ray fit is consistent with the observed one. This, though, is 
not strictly required for the consistency of the approach, since 
the CRE spectrum is not used to determine the diffusion pa- 
rameters. 

We also note that the nucleon injection spectra and the dif- 
fusion parameters are solely determined by fitting the local 
CR density, and neglecting the effect of DM. Injection of large 
quantities of nucleons from DM near the Galactic Center, that 
could alter the local abundances of CR nuclei (like the proton 
spectrum and the local B/C ratio, used to perform the fit) are 
in fact strongly excluded by anti-matter measurements. The 
anti-proton fraction, for example, rece ntly measured by the 
PAMELA experiment up to --100 GeV ( jAdriani et al.|2009b| l 
is about 10"^ above 10 GeV and thus constrains the DM con- 
tribution to be below this value. 

A drawback of using CR data is that they are usually af- 
fected by large systematics (for example, in the energy scale 
or in the solar modulation correction), so that the errors on 
the inferred diffusion parameters ar e larger than the st atisti- 
cal errors obtained from the fit (see [Trotta et al.[ ( [201 1| for a 
recent attempt to take into account these systematic effects). 
We checked that even quite large variations in these parame- 
ters affect our results only weakly. A more detailed discussion 
is defeiTed to Sec. [9] 

2.1. Limitations of the model 

The model described above represents well the gamma-ray 
sky, although various residuals (at a ^ 30% level ( [Ackermami 



et al.|201^) )), both at small and large scales, remain. ITiese 



residuals can be ascribed to various limitations of the mod- 
els: i) imperfections in the modeling of gas and ISRF com- 
ponents, ii) simplified assumptions in the propagation set-up 
(e.g. assumption of isotropy and homogeneity of the diffusion 
coefficient), iii) unresolved point sources, which are expected 
to con tribute to the diffuse emission at a level of 10% (Strong] 
2007 1, and have not been ta ken into account in the modeling 
m Ackermann et al. ( 2012a 1, iv) missing structures like Loop 



More precisely, the zi,=2, 15 cases are not reported in [Ackermann et ah] 
j20T2a|. The values used here for zh = 2, 15 are, however, obtained in the 
same way as the other z;, cases. 
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Parameter 


Value 


Halo Height z/,(kpc) 


2 


4 


6 


8 


10 


15 


Diffusion Coefficient Dq (cm^s"') 


2.7 X lO^** 


5.3 X lO^s 


7.1 X 10-** 


8.3 X 10^** 


9.4 X 102« 


1.0 X lO^' 


Diffusion Index 5 


0.33 


0.33 


0.33 


0.33 


0.33 


0.33 


Alfven Velocity (km s"') 


35.0 


33.5 


31.1 


29.5 


28.6 


26.3 


Nucleon Injection Index (Low) 7,, 1 


1.86 


1.88 


1.90 


1.92 


1.94 


1.96 


Nucleon Injection Index (High) 7,, 2 


2.39 


2.39 


2.39 


2.39 


2.39 


2.39 


Nucleon break rigidity pi,r,p(GV) 


11.5 


11.5 


11.5 


11.5 


11.5 


11.5 



Table 1 

CR diffusion parameters from |Ackermann et al.H2012a^ used in this work. 



I ("Casandji an et al.|2009 1 or the Galactic Bubbles/Lobes ( Su 
pal. 2010, ). 

We will deal with these limitations in various ways. Struc- 
tures like Loop I and the Galactic Bubbles appear mainly at 
high Galactic latitudes and their effects on the fitting can be 
limited using a Region of Interest (ROI) with limited extent in 
Galactic latitude. In the following we will consider a ROI in 
Galactic latitude, b, of 5° < \b\ < 15°, and Galactic longitude, 
Z, |Z| < 80°, see Sec. |5] As for small scale residuals, we be- 
lieve they are due to imperfections in the gas maps. In order to 
quantify their effect on the constraints of the DM properties, 
we will calculate the likelihood for several different assump- 
tions on the gas total column density. 

Despite the various choices described above and the large 
freedom we leave in the model (CRSDs, z/„ index of the elec- 
tron injection spectrum), we still see residuals in our ROI at 
the ±30% level and at > 3 cr significance (see the figures 
in Sec[8]l. Positive residuals, in particular, appear in vari- 
ous places in the ROI, especially in connection with the low 
Galactic latitude extension of the Lobes and Loop I. Residuals 
at the same level appear also when we perform a fit in a high 
latitude (I/?! > 40°) control region, where DM is not expected 
to contribute significantly, and thus seem to indicate the gen- 
eral level of accuracy achievable with the present modeling 
of the diffuse emission. Note that the residuals related to the 
Lobes and Loop I do not appear in the official Ferm/-LAT dif- 
fuse modej^ since there they are explicitly modeled through 
the use of patches. Since the residuals do not seem obviously 
related to DM, as it would be, for example, in the case of a sin- 
gle strong positive residual near the Galactic Center, we thus 
decide to focus in the following on setting limits on the possi- 
ble DM signal, rather than searching for a DM signal. Given 
the presence of the residuals like those described above, we 
also decide to quote more generous limits at the 3 and 5 sigma 
level, as well as conservative limits without assumptions over 
the astrophysical background. A dedicated search for a DM 
signal will be reported in a forthcoming paper, where data in 
the Galactic Plane and at low energies ( < 1 GeV) also will 
be exploited in order to help separating a real DM signal from 
other astrophysical processes which may be responsible for 
the above residuals. 

3. OUTLINE OF THE APPROACH TO SET DM LIMITS 

Having a parameterized model of the astrophysical Galac- 
tic diffuse emission as described above, we could imagine ex- 
ploring potential DM contributions by adding the DM com- 
ponent in the fit and performing a global joint fit of DM and 
the astrophysical model parameters. In practice, however, it is 
at the moment computationally challenging to perform such a 



large global fit. Thus, some additional simplifying assump- 
tion must be made. As akeady outlined in the previous sec- 
tion, the m ain simplification which we will adopt (which is 
also used in |Ackermann et alT ( |2012a |l) is the splitting of the fit 
in two separate steps constraining some parameters based on 
measurements of the local CR intensities while constraining 
other parameters based on Fermi LAT data. For the first step 
we rely entirely on Ackermann et al. (2012a) and we use the 
6 different diffusion models (for 6 different z/,) whose param- 
eters are summarized in Table [T] We then use only gamma- 
ray data to perform the fit for the CRSDs, zh, the electron 
index, gas maps with different column densities of the inter- 
stellar gas, and DM. Some additional parameters which enter 



the gamma-ray fit are further introduced in Sec. 7.3 

Our aim is to constrain the DM properties and treat the pa- 
rameters of the astrophysical diffuse gamma-ray background 
as nuisance parameters. Those parameters are typically cor- 
related with the assumed DM content and it is thus important 
to consider them since they affect directly the DM fit. It is 
clear, for example, that the CRSD should have a large influ- 
ence on the fit of the DM component. This can be seen from 
Figure [T] which shows that the gamma-ray signal produced 
by DM is somewhat degenerate with the IC signal from CR 
sources placed in the inner Galaxy. Besides small morpho- 
logical differences they mainly differ in the energy spectrum, 
which, however, is quite model dependent in the DM case. To 
explore the effect of uncertainties in the foreground model- 
ing we use the LAT data to fit the CRSD for both the nuclei 
and CREs, as well as the CRE injection spectrum which di- 
rectly affects the IC component. It is important to stress that 
to constrain DM in a self-consistent way, we will fit the CR 
distributions and DM at the same time, in order to take into 
account the degeneracy between the two. Fortunately, in our 
ROI, above and below the inner Galaxy, and around few GeV 
in energy, the tt*^ component, which has different morphol- 
ogy and is not degenerate with the DM signal, is dominant 
over IC and bremsstrahlung by about a factor 4-5. So, we 
expect approximately the same factor in improvement in DM 
constraints with respect to the fits that obtain limits without 
modeling this background. We will see in Sec. |6] and [8] that 
this expectation approximately holds. 

With the general approach above we set D M li mits using the 
profile likelihood method (outlined in Sec. |7.1| i. Besides the 
approach above, we will also quote conservative upper limits 
using the data only (i.e. without performing any modeling 
of the astrophysical background) . These conservativ e limits 
are along the lines of the work of Papucci & Strumia (2010); 
ICirelh et al.| ( |20I()] ), which use a similar approach to set DM 
limits based on the first year of Fermi LAT data. 



' For a description see http://fermi.gsfc.nasa.gov/ssc/data/access/lat/BackgroundModels.html. 



4. MAPS 
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4.1. DM maps 

The template maps used in the fits to model the DM con- 
tribution depend on the assumed DM distribution and the as- 
sumed annihilation/decay channel. Numerical simulations of 
the Milky Way-scale halos indicate a smooth distribution that 



contains a large number of subhalos (Diemand et al. 



2007b 



|Springel"et al. 2008b ). The gamma-ray signal from the sub 
halo population is expected to dominate in the region of the 
outer halo, while in the inner <20° region of the Galaxy, its 
contribution is expected to be subdominant (Diemand et al. 
2007al[Sprmgel et al.|2008aHPieri et al.l20Tr) ). In our ROl the 
subhalo contribution therefore should be mild and we conser- 



vatively consider only the smooth component in this work. 
We parametrize the smooth DM density p with a NFW spatial 
profile (Navarro et al. 1996 ) and a cored (isothermal-sphere) 
profile ( ,Begeman et al.,1991,,Bahcall & Soneira,1980| ): 



NFW: p(r) = po 1 + 



Re 



1 



Isothermal : p(r) = po 



(1) 



(2) 



These are traditional benchmark choices, as NFW is moti- 
vated by N-body simulations, while cored profiles are instead 
motivated by the observations of rotation curves of galax- 
ies and are also found in simulations of a Milky Way-scale 
halos involving baryons (Maccio et al. 20l2). The Einasto 
profile ( Merritt et al. 2006, Navarro et al._2 010) is emerging 
as a better fit to more recent numerical simulations, but for 
brevity we do not consider it here. It is expected that this 
profile should lead to DM limits stronger by ^ 30% in our 
ROI, with respect to a choice of a NFW profile ( [Cirelli et al.| 
|2010| l. The main uncertainty in the DM halo profile comes 
from the poorly known (and modeled) baryonic effects. In- 
deed, with our choice of NFW and Isothermal profiles, our 
aim is to roughly bracket the uncertainties expected from the 
DM profile. For the lo cal density of DM w e take the value 
of po = 0.43 GeV cm'^ (jSalucci et al.|2010| and the scale 
radius is assumed to be Rs = 20 kpc (NFW) and Rc = 2.8 kpc 
(isothermal profile). The actual choice of the DM density pro- 
file does not have a major effect on our limits (see Sec. [8]l 
as we do not consider the central few degrees of the Galaxy 
(where these distributions differ the most). A choice of a more 
extended core of ^ 5 kpc seems possible, altho ugh less fa- 
vor ed by data (Bergstrom et al . 1998) (s ee also, Catena & Ul-| 
[Iiol ( poTUl l; [Weber & de Boer| (2010); l occo et al.| ( |2011| l for 
further discussions on the po ™d DM profile uncertainties). 
With this choice our limits would worsen by a factor of < 2. 
We also set the distance of the solar system fro m the center 
of the Galaxy to the value Rq = 8.5 kpc (Kerr & Lynden-Bell 
[1986) 1. 

For the annihilation/decay spectra we consider three chan- 
nels with distinctly different signatures: annihilation/decay 
into the bb channel, into p.'^pT, and into t^t~ . In the first case 
gamma rays are produced through hadronization and pion de- 
cay. The resulting spectra are similar for all channels in which 
DM annihilations/decays produce heavy quarks and gauge 
bosons in the energy range considered here and is therefore 

" The measurement has a typical associated error bar of ifcO. 1 GeV cm~^ 
and a possible spread up to 0.2-0.7 GeV cm" dSalucci et al. 2010 CireUi 
|etal.|201l| . ^ ' " ' 



representative for a large set of WIMP particle physics mod- 
els. The choice of annihilation/decay into leptonic channels, 
provided by the second and third scenarios, is motivated by 
the PAMELA positron fraction ( Adriani et al. 2009a) and the 
Fermi LAT electrons plus positrons (Abdo et al. 2009 ) mea- 
surements (for interpre tation of thes e measur ements in terms 
of a DM signal see e.g. [Grasso etal 2009 ; M eade et al.|2010 



Bergstrom et al.[|2009 |). In this case, gamma rays are domi 



nantly produced through radiative processes of electrons, as 
well as through the Final State Radiation (FSR). 

We produce the DM maps with a version of GALPROP 
slightly modified to implement custom DM profiles and in- 
jection spectra. For the prompt photons case GALPROP inte- 
grates along the fine of sight the DM gamma-ray emissivity 
given by Q^{r,E) = p^ {(Tav) /2?m^ x dN^/dE in the annihila- 
tion case and by Q^(r,E) = pTo/m^ x dN^/dE in the decay 
case, m^ being the DM particle mass, dN^/dE the gamma- 
ray annihilation/decay spectrum, {(Jav) the thermally aver- 
aged DM annihilation cross section and = l/r, the DM 
decay rate (i.e. the inverse of the DM lifetime). In the cases 
where the propagation of electrons produced in DM annihi- 
lation is relevant to the resulting gamma-ray emission, tem- 
plate maps are produced including the IC emission from DM 
annihilation/decay-generated CREs, which have been propa- 
gated using the relevant set of propagation parameters. An 
example of the resulting DM maps at 10 GeV for two DM 
models (annihilation to bb and p^p~ channels) is shown in 
Figure [1] for a DM mass of 250 GeV. 

We calculate the DM injection spectrum of electrons and 
gamma rays by u sing the PPPC4DMID tool described in 
[Cirelh et"ar] ( |201 1| . This package provides interpolating func- 
tions calculated from a si mulation o f DM annihilation/decay 
with the PYTHIA (Sj ostrand et al.[r2008 ) event generator. It 
has recently been shown that for DM candidates with masses 
above the electroweak scale, the standard predictions for the 
annihilation/decay spectra are altered when account is taken 
of the production of electroweak gauge bosons from the FSR 
(&icheMe6 et al. 2009; Ciafaloni et al. 2011b a). We include 



these electro-weak corrections here, following Cirelli et al. 
( [201 1[ ), even though the effect on DM limits is marginal for 
our choice of energy range. 

4.2. Gamma-ray emission from CR interactions with 
interstellar gas 

The astrophysical diffuse emission predicted by GALPROP 
depends on assumptions about the distribution and column 
density of the interstellar gas entering the model. A signifi- 
cant uncertainty is related to the total gas column density due 
to the presence of dark gas. The dark gas contribut ion is es- 
timated using the Schlegel, Finkbeiner, Davis (SFD) [Schlegel| 

fal. ( 1998 ) E(B-V) dust reddening map as a tracer of the to- 
integrated gas column density (dark -i- radio-visible). The 
intensity of the velocity resolved H I maps derived from the 
21 -cm survey data is then rescaled by the ratio between the 
total and radio-visible column densities to account for dark 
gas. This rescaling factor depends on the ratio between the 
dust and gas column densities and thus on the dust to H I ratio 
(d2H I) and the dust to CO ratio (d2CO). In Ackermann et al.| 
(2012a) the dust to H I and CO gas ratios are fixed through a 
regression procedure of the SFD map to the radio-derived gas 
maps (see [Ackermann et al.[ ( |2012a| ) for more details) which 
yields d2HI=0.0137 x lO"^" mag cm^ and d2CO=0.0458 mag 
(K km s"')"'. Here we will instead explore different values 
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of the dust to gas ratios, to consider the possibility that this 
quantity is diffe rent in our ROI with resp ect to the all-sky de- 
rived value of Ackermann et al. (2012a) and thus, possibly, 
improv e the residuals, especially at small scale, as discussed 
in Sec. |2.1| In particular we will use 6 values of d2HI equally 



spaced in the range (0.0120-0.0170) xlO"^" mag cm^. The 
d2CO ratio will be instead fixed to the value 0.04 mag (K km 
s"')"'. This is justified in the light of the fact that above ±5° 
of Galactic latitude, where we perform the fit, there is very 
little CO and we are thus not very sensitive to this parameter. 
Nevertheless, we checked the results for different values of 
d2C0 in the range 0.03 - 0.06 mag (K km s"' )"' and found no 
appreciable change in the results. 

The E(B-V) map is not a reliable tracer of total column den- 
sity when the dust reddening becomes very high. Thus a cut in 
E(B-V) needs to be employed to exclude the regions of high 
Galactic extinction. In Ackermann et al. (2012a) a value of 
E(B-V)< 5 mag is found to be adequate and we will use it in 
the following. This cut excludes from the dark gas correction 
a narrow strip of few degrees along the Galactic Plane, which 

is a region outside our ROI. 

The spin temperature of the H I, Ts (see | Ackermann et al. 



( |2012a ) for the detailed definition), used to extract the H I 
maps from the radio data is also not very well kn own and in- 
troduces further uncertainties. However, since inl Ackermann| 



et al. ( 2012a i the total column density of the ISM is estimated 
from dust, the effect of is typically subdominant with re- 
spect to the dark gas correction. We will thus not consider it 



and fix 7^ to the typical value Ts = 150 K (Ackermann et al. 
|20T2a| . 

5. DATA SELECTION AND REGION OF INTEREST 

We use 24 months LAT data starting from 2008 August 5 
to 2010 July 31, in the energy range between 1 GeV and 100 
GeV. However, we use energies up to 400 GeV when deriv- 
ing DM limits with no assumption on the astrophysical back- 
ground, see next section. As the data above 100 GeV have 
poor statistics they have little weight when performing the fit 
to the gamma-ray data to determine the background model. 
Instead, for the no-background limits, where no fitting is per- 
formed, high energy points are relevant to determine the DM 
limits only for the hard FSR spectrum. We use only events 
classified as gamma rays in the P7CLEAN event selection and 
the corresponding P7CLEAN_V6 instrument response func- 
tion^ The data have bee n extracted and processed w ith the 
Fermi tools as described in [Ackermann et aT1 ( |2012a) . In or- 
der to minimize the contribution from the very bright Earth 
limb, we apply a maximum zenith angle cut of 90° . In addi- 
tion, we also limit our data set to include only photons with 
an incidence angle from the instrument z-axis of < 72°. The 
events are divided into 5 logarithmically spaced energy bins. 
The total number of photons in our ROI is ^ 350000. We 
use a HEALPiji ^ ' ( Gorski et al. 2005 1 n s i de=64 pixelization 
scheme for the spatial binning, corresponding to a bin size of 
approximately 0.9 deg x 0.9 deg. We f urther mask the po int 
sources from the IFGL source catalog ( Abdo et al.l2010bll to 



limit the impact of bright point sources on the fit. We mask all 
pixels within 1 deg of the location of the point sources. With 
this choice we remove ^ 25% of the photo ns, which leaves us 
with ~ 270000. The 2FGL source catalog ( [Nolan et al.|20T2] l 
contains some further weak point sources which we do not 

http://fermi.gsfc.iiasa.gov/ssc/ 
^' http://healpix.jpl.nasa.gov 



mask to not impact dramatically the size of our ROI. The ef- 
fect of masking 2FGL sources is further considered in Sec. [8] 
Since we focus on the Galactic halo we choose a ROI limited 
by ±15° in Galactic latitude and ±80° in longitude, to focus 
on the region wh ere the S/N ratio for DM is the highest ( |Ser^ 
pico & Zaharijas 2008), and to minimize the effect of the nigh 
latitude structures like Loop I or the Galactic Lobes. This ROI 
also excludes the Outer Galaxy where the DM signal should 
be lower. Furthermore, we mask the region \b\ <5 deg along 
the Galactic Plane, in order to reduce the uncertainty due to 
the modeling of the astrophysical and DM emission profiles 
discussed above. 

6. DM LIMITS WITH NO ASSUMPTION ON THE 
ASTROPHYSICAL BACKGROUND 

To set DM limits with no assumption on the astrophysi- 
cal background we first convolve a given DM model with 
the exposure maps and Point-Spread Function (PSF) to obtain 
the counts expected from DM annihilation. The GaRDiAn 



(Gamma-Ray Diffuse AN alysis) software package (Acker- 
[mann et al.| 2012a[ [2009 | l is used for this processing^ The 
expected counts are then compared with the observed counts 
in our ROI and the upper limit is set to the minimum DM 
normalization which gives counts in excess of the observed 
ones in at least one bin. More precisely, we set 3 a upper Um- 
its giv en by the requirement n,x)M-3^n,x)M > n,- dCowan et al. 
[201 l| l, where n,oM is the expected number of counts trom UM 
in the bin / and n, the actual observed number of counts. It 
should be noted that the formula assumes a Gaussian model 
for the fluctuations, which is a good approximation given the 
large bin size and the number of counts per bin we use in 
this case (see below). The large Poisson noise present espe- 
cially at high (> 10 GeV) energies due to the limited number 
of counts per pixel, affects the limits for DM masses above 
100 GeV, weakening them somewhat. To reduce the Pois- 
son noise, only for the present case of no background mod- 
eling we choose a larger pixel size so to increase the number 
of counts per pixel. However, a very large pixel size would 
wash out the DM signal, diluting it in large regions, again 
weakening the limits. We chose the case with a pixel size of 
about 7° X 7° (nside=8) since it gives a reasonable com- 
promise between the two competing factor^ In this way 
limits typically improve by a factor of a few with respect to 
the case nside=64. Limits for DM masses below 100 GeV, 
instead, are only very weakly affected by the choice of the 
pixel size in the range l°-7°. Finally, again only for the 
present case of no modeling of the background, we use an 
extended energy range up to 400 GeV. This, in practice, is im- 
portant only for the case for masses above 100 GeV and 
when we consider FSR only (since the ji^jT FSR annihilation 
spectrum is peaked near the energy corresponding to the DM 
mass and thus can be constrained only by using higher-energy 
data). For the other cases, instead, there is always significant 
gamma-ray emission below 100 GeV, either from prompt or 
IC photons and the extended energy range does not affect the 
limits appreciably. 

The limits derived from this analysis are discussed in SecjS] 
These constraints are about a factor of 5 worse than those ob- 
tained with a modeling of the background (see next section), 
which is in agreement with the estimate made in Sec. [3] 

The mask is always defined (and applied) at n s i de = 64. After applying 
the mask the data (and the models) are downgraded to the larger pixel size. 
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Figure 2. Example of profile likelihood curves for four different DM annihilation/decay scenarios. Each curve refers to a particular model of the background. The 
envelope of the various curves approximates the global profile likelihood marginalized over the astrophysical uncertainties accounted for in our fitting procedure. 
The curve corresponding to the model setting the global minimum, y„„„, is highlighted in red. The y scale is arbitrarily re-shifted so that the minimum value is 
zero. The green curve corresponds to the model setting the 3 a upper limit (i.e. the model which is both part of the envelope profile likelihood and intersects the 
horizontal line located at +9). The upper limit is then effectively given by the x coordinate of the intersection point. The blue curve is similar, but for the 5 a case 
(and intersects the horizontal line located at -F25). For these 3 models the corresponding values of z/,, 7,, 2, and d2HI are given in the caption. Panel description: 
10 GeV DM particle decaying (DEC) into bb and NFW profile (upper left), 91 GeV DM particle annihilating (AN) into bh and NFW profile (upper right), 5 GeV 
DM particle decaying into t+t" and NFW profile (lower left) and 750 GeV DM particle annihilating into r+r" and NFW profile (lower right). 



7. DM LIMITS WITH MODELING OF 
ASTROPHYSICAL BACKGROUND 

We derive a second set of upper limits taking into account 
a model of the astrophysical background. As described in 
sections[2]and[3] the approach we use is a combined fit of DM 
and of a parameterized background model and we consider 
the uncertainties in the background model parameters through 
the profile likelihood method described below. 

7.1. Profile Likelihood and grid scanning 

For each DM channel and mass the model which describes 
the LAT data best maximizes the likelihood function which is 
defined as a product running over all spatial and spectral bins 
i , 

LkiOoM) = LkiOoM, a) = maxsY[Pik(ni;d, Oom), 0) 

i 

where Pik is the Poisson distribution for observing n, events 
in bin / given an expectation value that depends on the pa- 
rameter set {9dm, 0). 9dm is the intensity of the DM compo- 
nent, a represents the set of parameters which enter the as- 
trophysical diffuse emission model as linear pre-factors to the 
individual model components (cf. equation [5j below), while 
k denotes the set of parameters which enter in a non-linear 
way. Individual GALPROP models have been calculated for 
a grid of values in the k parameter space. For each family of 



models with the same set of non-linear parameters k the pro- 
file likelihood curve is defined for each 9dm as the likelihood 
which is maximal over t he possible choices of the parame- 
ters a for fixed 9dm (see Rolke et al. ( 2005| l and references 

therein). The notation a represents the conditional maximiza- 
tion of the likelihood with respect to these parameters. The 
linear part of the fit is performed with GaRDiAn, which for 
each fixed value of 9dm finds the a parameters which maxi- 
mizes the likelihood and the value of the likelihood itself at 
the ma ximu rrp^ (for details about fitting linear parameters see 
7.3 I. However, since building the profile likelihood 



section 

on a grid of 9dm values is computationally expensive, we use 
an alternative approach including 9dm explicitly in the set of 
parameters fitted by GaRDiAn. In this case GaRDiAn also 
computes the 9dm value which maximize the likelihood (the 
best fit value 9dmo) and its la error estimated from the cur- 
vature of the logL^^ around the minimum. We then approxi- 
mate the profile likelihood as a Gaussian in 9dm with mean 
9dmi) and width ctsomo • We have verified that this approxima- 
tion works extremely well for a subset of cases for which we 
also explicitly computed the profile likelihood, tabulating it 
on a grid of 9dm values. We will thus use this approximation 

'-^ Technically, instead of maximizing the likelihood, GaRDiAn minimizes 
the (negative ot^ log-likelihood, -logL, using an external minimizer. For our 
analyses we used GaRDiAn with the Minuit l |James & Roos|1975f mini- 
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throughout the rest of the analysis. 

In this way we end up with a set of k profiles of likelihood 
LkifioM), one for each combination of the non-linear param- 
eters. The envelope of these curves then approximates the 
final profile likelihood curve, L{Oom), where all the parame- 
ters, linear and non-linear have been included in the profilq^ 
. Examples of such final profile likelihood curves for specific 
DM models can be seen in Figure |2j and will be discussed 
more in detail in Sec. 17.31 

Limits are calculated from the profile likelihood function by 
finding the 9r,M.iim values for which LieoMMm)/UdoMjna.x) is 
exp(-9/2) and exp(-25/2), for 3 and 5 a C.L. limits, respec- 
tively. This approximation is exact for Gaussian likelihood 
functions in one parameter and, due to invariance of the like- 
lihood function under reparameterization, it is most often also 
applicable to the non-Gaussian case (James 2006). For the 
case of handling nuisance parameters, this is not true a pri- 
ori, but has been shown to give satisfactory propertie s for a 
variety of nuisance parameter configurations (e.g. in |Rolke 



[etiL] ( |2005l l; |Abdo et al.|P010c| ); |Ackermann et al.| ( |2011| i) 
In particular see also the recent search for the Higgs boson at 
the Large Hadron Collider, where 0(100) nuisance p arame- 
ters need to be taken into account (ATLAS Collaboratio n Aad] 
et al. 2012|. We therefore are confident that this approach 



gives the desired statistical properties, i.e., good coverage and 
discovery power, also in our analysis. 

7.2. Free CR Source Distribution and constrained setup 

limits 

In this section we introduce the first set of linear parameter, 
i.e. the coefficients defining the CRSDs. The remaining linear 
parameters will be introduced in the next section. 

As noted in section |2l C RSDs ( for example the ones con- 
sidered in Acker mann et ar] ([ 20 12a )) can be modeled from the 
direct observation of tracers of SNR, and so can be observa- 
tionally biased. The uncertainty in the distribution of the trac- 
ers in the inner Galaxy is therefore large and should be taken 
into account in the derivation of the DM limits. We therefore 
fit the CRSD from the gamma-ray data, as described below. 

Due to the linearity of the propagation equation it is pos- 
sible to combine solutions obtained from different CRSDs. 
To exploit this feature we define a parametric CRSD as sum 
of step functions in Galactocentric radius R, with each step 
spanning a disjoint range in R: 



e,pCRSD(R) = 0(R-Ri)0{Ri+i -R) 



(4) 



We choose 7 steps with boundaries: Rj =0, 1.0, 3.0, 5.0, 7.0, 
9.0, 12.0, 20.0 kpc. The expected gamma-ray all sky emission 
for each of the 14 single-step primary e and p distributions are 
calculated with GALPROP. It is also worth noting that a dif- 
ferent GALPROP run needs to be done for each set of values 
of the non-linear parameters, since, for a given e,/:iCRSD the 
output depends on the entire propagation setup. For more ac- 
curate output, especially in the inner Galaxy, which we are 
interested in, GALPROP is run with a finer grid in Galacto- 
centric radius R with dr = 0.1 kpc, compared to the standard 
grid of dr=l kpc. The coefficients cf are set to unity for the 
individual GALPROP runs and then fitted from the gamma-ray 
data as described below. 

We will sometime use in the following the term marginalizing although, 
typically, the term applies only within the framework of Bayesian analyses. 
In our frequentist approach it is called profiling. 



In order to have conservative and robust limits we constrain 
the parameter space defined above by setting c\'' = c^' = 0, 
i.e. setting to zero the e,/?CRSDs in the inner Galaxy region, 
within 3 kpc of the Galactic Center In this way, potential 
e and p CR sources which would be required in the inner 
Galaxy will be potentially compensated by DM, producing 
conservative constraints. A second important reason to set the 
inner e,/7CRSDs to zero is the fact that they are strongly de- 
generate with DM (especially the inner eCRSD, see Figure 
[l]). Besides slight morphological differences, an astrophysi- 
cal CRE source in the inner Galaxy is hardly distinguishable 
from a DM source, apart, perhaps, from differences in the 
energy spectrum. To break this degeneracy we would need 
to use data along the Galactic Plane (within ±5° in latitude) 
since these are expected to be the most constraining for the 
e,/?CRSDs in the inner Galaxy. However, the Galactic Center 
region is quite complex and modeling it is beyond the scope of 
the current paper. We therefore defer such a study to follow- 
up pubhcations. 

7.3. Fitting procedure 

In the fit of the expected gamma-ray emission to the Fermi 
LAT data we determine the normalizations of the contribu- 
tions from DM and from each step of the CRSD function de- 
fined in Eqn. |4]that best fit the data. To achieve this, we 
need to split each contribution into several components cor- 
responding to the type of target and physical process respon- 
sible for the emission. The emission from tt" decay depends 
only on the distribution of the CR nuclei sources, while the 
emission from bremsstrahlung and inverse Compton depends 
only on the distribution of the CR electron sources (emission 
from interactions of secondary electrons produced in CR nu- 
clei interactions is negligible above 1 GeV). The gamma-ray 
emission arising from interactions of CRs with molecular gas 
traced by CO depends further on the assumed conversion fac- 
tor Xco between the CO intensity and the column density of 
the molecular gas. This conversion factor is uncertain and we 
vary it freely for each annulus. We determine effective Xco 
factors implicitly in the fit by splitting the calculated expected 
gamma-ray emission from CR interactions with molecular gas 
into Galactocentric annuli which are separately normalized. 
Additionally an isotropic component arising from the extra- 
galactic gamma-ray background and misclassified charged 
particles needs to be included to fit the Fermi LAT data. We 
do not include sources in the fit as we use a mask to filter the 
IFGL point sources (cf. Sec [5]). To rule out the possibility 
that some bright sources mightleak out of the mask and bias 
the fit we performed test fits including explicitly the IFGL 
point sources as a further template map, finding that the inclu- 
sion of the point sources introduces only a negligible change 
in the results. Equation [5] summarizes how we parametrize 
the expected gamma-ray emission / in the fit based on the 
components mentioned above. Each component is calculated 
using GALPROP and is available as a template map after the 
GALPROP run. In summary, the various GALPROP outputs 
are combined as: 

j ^ 
(Xt + Xic) + E "/GB,m I GB" . (5) 



11 



Non linear Parameters 


Symbol 


Grid values 


index of the injection CRE spectrum 
lialf lieiglit of the diffusive halcj^ 
dust to HI ratio 


7c.2 
d2HI 


1.925, 2.050, 2.175, 2.300, 2.425, 2.550, 2.675, 2.800 
2, 4, 6, 8, 10, 15 kpc 
(0.012, 0.013, 0.014, 0.015, 0.016, 0.017) X lO^^o mag cm^ 


Linear Parameters 


Symbol 


Range of variation 


eCRSD and pCRSD coefficients 
local H2 to CO factor 
IGB normalization in various energy bins 
DM normalization 


yloc 
aiGB.m 


0,+oo 

0-50 X lO^o cm-- (K km s^'r' 
free 
free 



"The parameters Dq, 5, va, 7p.i, 7p.i, Phr.p are varied together with zh as indicated in Table[T] 



Table 2 

Summary table of the parameters varied in the fit. The top part of the table shows the non linear parameters and the grid values at which the likelihood is 
computed. The bottom part shows the linear parameters and the range of variation allowed in the fit. The coefficients of the CRSDs are forced to be positive. 



except c^"' and c^'' which are set to zero. The local Xco ratio is restricted to vary in the range 0-50 X 10^" cm - (K km s 

to assume both positive and negative values. See the text for more details. 



) , while aioB.m and are left free 



\ w/o DM 




\ bb, AN 




\ mu, AN 




\ tau, AN 




\ \ bb, DEC 




, \ \ rau, DEC 




\ \ \ tim. DEC 




**X \ \ NEW 




V \ ISO 









- w/o DM 

- bb, AN. NEW 

- mil. AN. NEW 

- tau, AN, NEW 
■ bb, DEC, NEW 

mu, DEC, NEW 
tau, DEC, NEW 



2.2 2.3 
ye.2 




0.013 0.014 0.015 0.016 

d2HI[10-='° magcm'^] 



- w/o DM 

- bb, AN, NEW 

- mu, AN, NEW 

- lau, AN, NEW 



- bb, DEC, NEW 
■ mu, DEC, NEW 

- tau, DEC, NEW 




8 

z [kpc] 



Figure 3. Profile likelihood curves for z;,, 7^ 2 and d2HI. The various curves refer to the case of no DM or different DM models (see the legend in the figure, 
where we mark a dominant decaWDEC) or annihilation (AN) channel and the assumed DM profile). All minima are normalized to the same level. Horizontal 
dotted lines indicate, as in Figure|2] a difference in — 2AlogL from the minimum of 9 (3a) and 25 (5cr). 



The sum over / is the sum over all step-like CRSD functions, 
the sum over j corresponds to the sum over all Galactocen- 
tric annuli (details of the procedure of a placement of the 
gas in Galactocen t ric annu li and their boundaries are given in 
Ackermann et al. ( 2012a 1). H denotes the gamma-ray emis- 
sion from atomic and ionized interstellar gas while H2 the one 
from molecular hydrogen and IC the Inverse Compton emis- 
sion. X7 and Xic are the prompt and Inverse Compton (when 
present) DM contribution and the overall DM normal- 
ization. IGB"' denote the Isotropic Gamma-ray Background 
(IGB) intensity for each of the five energy bins over which 



the index m runs. For better stability of the fit the template for 
IGB'" is build starting from an 1GB with a power law s pec- 
trum and normalization as given in Abdo et al. ( 2010e| l. In 
this way the fit coefficients ajcB.m are typically of order 1 . In 
all the rest of the expression in Eqn. ISlthe energy index m is 
implicit since we don't allow for the freedom of varying the 
GALPROP output from energy bin to energy bin. Finally, it 
should be also noted that in our case, where we mask ±5° 
along the plane, the above expression actually simplifies con- 
siderably since only the local ring Xco factor enters the sum, 
since all the other II2 rings do not extend further than 5 de- 
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grees from the plane. Also to be noted is the fact that, since in 
Eqn. |5]//2 denotes a gamma-ray emission map, the expression 
has been already intrinsically multiplied by an Xco factor to 
convert the CO line intensity into an H2 column density. We 
in fact normalize all the H2 gamma-ray maps using the value 
Xco = 1 X 10-° cm-2 (K km s"')"'. The Xco in Eqn. |5]are 
thus adimentional ratios with respect to the reference value 
1 X 10^° cm"^ (K km s"')"'. With a slight abuse of notation 
we denote them also as Xco factors. 

The above expression predicts the expected gamma-ray 
counts in terms of the parameters (cf', X^q, ajcB.m and 
for a total of 7h-7h-1h-5h-1=21 parameters). GaRDiAn is used 
to build the profile likelihood for the intensity of the DM com- 
ponent by finding the set of parameter values which max- 
imize the likelihood for a given . 

The outlined procedure is then repeated for each set of val- 
ues of the non-linear propagation and injection parameters to 
obtain the full set of profile likelihood curves. We scan over 
the following three parameters: the half-height of the diffu- 
sive zone Zh, the index of the electron injection spectrum j^.i 
and the dust-to-H I ratio d2HI. Specifically, we choose 6 val- 
ues of z/, = 2,4,6,8, 10, 15 kpc, 8 values of -fe,! linearly spaced 
between 1.925 and 2.8, and 6 values of d2Hl linearly spaced 
in the range (0.0120 - 0.0170) x 10"^° mag cm^. Taking into 
account the 7 step functions used for the e,pCRSDs we scan 
over a grid of 7x5x8x6=1680 GALPROP models (or rather 
GALPROP runs since combinations of the steps are effectively 
a single GALPROP model.). 

In order to follow more easily the entire fitting procedure 
we report in Table [2] a summary of all the parameters em- 
ployed in our analysis, linear and non-linear, together with 
their range of variation in the fit or discrete values used in the 
grid. 

Figure|2]shows some examples of the profile likelihoods for 
selected DM masses and annihilation channels. The limits are 
set by first finding the absolute minimum and then looking at 
the intersection between the envelope of the various parabolae 
and the 3 and 5 cr horizontal lines. An important point to note 
is that, for each DM model, the global minimum we found lies 
within the 3(5) a regions of many different models. This is a 
basic sanity check against a bias in our procedure, as would be 
suspected if the model giving the minimum was inconsistent 
with the bulk of the other models considered. This point is 
further illustrated in Figure|3] where the profile likelihoods for 
the three nonlinear parameters, z/,, 7f 2 and d2Hl, are shown. 
To ease reading of the figure the profiling is actually per- 
formed with further grouping DM models with different DM 
masses, but keeping the different DM channels, DM profiles 
and the annihilation/decay cases separately. The curve for the 
fit without DM is also shown for comparison. Each resulting 
curve has been further rescaled to a common minimum, since 
we are interested in showing that several models are within 
-2AlogL < 25 around the minimum for each DM fit. The 7^ 2 
profile, for example, indicates that all models with 7^, 2 from 
1.9 to 2.4 are within -2AlogL < 25 around the minimum il- 
lustrating that the sampling around each of the minima for 
the six DM models is dense. Similarly, the d2Hl profile indi- 
cates that all models with d2HI in the range (0.120 - 0.160) 
X 10"^" mag cm^ are within Sa from the minima for each of 
the six DM models. Finally the zi, profile indicates that basi- 
cally all the considered values of zh are close to the absolute 
minima. This last result is not surprising since, within our 
low-latitude ROI, we have little sensitivity to different Zh and 



basically all of them fit equally well. There is some tendency 
to favor higher values of Zh when DM is not included in the 
fit, while with DM the trend is inverted, although the feature 
is not extremely significant it is potentially very interesting. 

As explained in Sec. [3] in our analysis the DM parameter is 
the one of prime interest and we thus treat the parameters (lin- 
ear and non-linear) of the diffuse emission as nuisance param- 
eters which we include to take into account degeneracies with 
DM (i.e to marginalize over them) and establish more robust 
limits. This is a reasonable assumption since, for example, 
most of the IC emission and thus the sensitivity to 7^, 2 comes 
from data within 5° of the Galactic Plane, which are not con- 
sidered here. Similarly, to be sensitive to z/,, high Galactic lat- 
itude data should be included. These plots thus should be re- 
garded only as indicative of the achievable constraints, while 
a careful analysis will be deferred to later publications. 

An issue that we are not addressing here is whether DM is 
required or not in the fit, and, in the former case, finding the 
best model among different DM models. Since we have seen 
that systematic uncertainties related to the limitations in mod- 
eling astrophysical contributions to the Galactic diffuse emis- 
sion are comparable in size to any DM signal we fit, we have 
focussed on setting constraints on potential DM contributions 
to the Galactic diffuse emission. The systematic uncertain- 
ties in the Galactic diffuse emission modeling likely could 
be reduced by including the Galactic Plane/Galactic Center 
data in the analysis. Furthermore, a realistic study of the 
problem would require also considering other possible com- 
ponents which might be present in the Halo, like contributions 
from a population of unresolved pulsars or the emission from 
the Bubbles/Lobes. We defer this analisys to a subsequent 
study. 

8. RESULTS 

Despite the various conservative choices described above, 
the resulting limits are quite stringent. Upper limits on the ve- 
locity averaged annihilation cross section into various chan- 
nels are shown in Figure |4] for an NEW and isothermal pro- 
file of the DM halo. The limits obtained without modeling of 
the astrophysical background are comp atible with the resul t 
of similar analyses presented in Papucci & Strumia (2010 1; 
CirelU et al. ( 2010[ l. Limits with model of the background, 
instead, are significantly improved with respect to the above 
ones. They are competitive with respect to the limits from 
LAT searches for a signal from DM a nnihilation/decay in 
dwarf galaxies (Ackermann et al.|2011| l and Galaxy clusters 
( |Ackermann et al.|2010| l. 

In particular, as shown in Figure |4] for masses around 20 
GeV the thermal relic value of the annihilation cross section 
is reached, both for the bb and t^t~ channels. The limits are 
also imp roved over the ones deri ved in the analysis of dwarf 
galaxies ( [Ackerm ann et al. 201 1) which did not consider the 
inverse Compton emission (in dwarfs this component is quite 
uncertain) and also improved over constraints imposed on DM 
annihilations from the absence of a measurable effect on CMB 
anisotropics ( [Galli et al.|2011| l. 

A limitation of our constraints is their dependence on 
poorly determined properties of the Milky Way DM Halo, 
in particular on po, from which the normalization of the DM 
signal, and thus the limits, depends quadratically in the an- 
nihilation case and linearly in the decay case. We use the 
recent determination po = 0.43 GeV cm"^ from Salucci et al. 
( |2010[ l, which has, however, a large uncertainty, with a typ- 
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ical associated error bar of ±0.1 GeV cm ^ 

-3 



and a possible 



Cirelli 



spread up to 0.2-0.7 GeV cm"^ ( [Salucci et al.|j2010[ , 
et al.|20lT i. Whether the limits will worsen, or improve, thus 
awaits a better determination of po- To show the effect of the 
Pq uncertainty on the limits we plot, for illustration, in the top 



left panel of Figs|4]|5]the uncertainty band (red dotted lines) in 
the 3(7 no-background limits which would result from varying 
the local DM density po in the range 0.2-0.7 GeV cm"^ (we 
conservatively take here the larger scatter to show the max- 
imal impact of the uncertainty of po). A similar band, not 
shown in the plot for clarity, would be present for the lim- 
its including a model of the astrophysical background. The 
band is likely a generous estimate of the uncertainty since the 
variation of po is typically correlated with other properties of 
the DM Halo, such as the density profile and the distance of 
the solar system from the Galactic center Rs- The uncertainty 
band shown should just be considered an illustration, while a 
detailed study would be required to address the actual uncer- 
tainty, which is beyond the scope of this work. 

In Figure |4| w e also sho w the regions of the parameter space 
derived in Cirelli et al.||2010) from a DM fit to the Fermi 
LAT electron/positron data and PAMELA positron fraction 
data. Contours ai-e shown at 95% and 99.999% CL. These 
regions are rescaled down of a factor (0.43/0.3)^ ^ 2 to take 
into account the different local DM density p^) used in the two 
works (for the same reason the regions in the decay case are 
rescaled up by a factor of (0.43/0.3) ^ 1.4). We must also 
take into account that different en ergy losses o f local CREs 
are used here when compared to |Cirelli et al. (2010j). For 
CREs with E >\ TeV diffusion can be approximately ne- 
glected so that, from the diffusion loss equation, the steady 
state CRE flux scales approximately linearly with the energy 
loss time scale r (see Bergstrom et al. (2009);|Borriello et al.| 
(j2009a|). The VAMEhAIFermi regions, instead, scale as the 
inverse of the local CRE flux. Beyond the different local mag- 
netic field and ISR F, a major difference with the analysis in 
ICirelU et aT] ( |2010[ ) is that they neglect Klein-Nishina atten- 
uation effects in the Inverse Compton cross section, which 
are instead taken into account in our GALPROP calculations. 
These effects make the energy losses of CREs with E>\ TeV 
almost negligible on the optical part of the ISRF, and slightly 
decrease the energy losses due the infrared part. CMB losses 
are instead unaffected. As result, the energy loss time scale, 
and thus the CRE flux, increases (the PAMELA/ferw/ regions 
will be pushed downward). On the other hand, we use a local 



magnetic field of 5 pG as opposed to the 3 ^G used injCirelli 
et al.| ( [2010) , and this has the contrary effect, decreasing the 



energy losses time scale (which scales as Br^) and thus push- 
ing up again the PAMELA/f erm/ confidence regions. Over- 
all, to correctly derive the positions of the PAMELA/ferm/ 
confidence regions, the uncertainties in the local ISRF and 
magnetic field should be taken into account and marginal- 
ized away, which we leave for a follow-up work. Here, from 
the above considerations, we estimate that the location of the 
PAMELA/Ferm/ regions has a further uncertainty of a fac- 
tor ~ 2, so that they can possibly touch the exclusion limits. 
Thus, we cannot robustly rule out the DM annihilation in- 
terpretation of the PAMELA/Ferm/ CRs anomalies, although 
this interpretation is challenged. Finally, we note that the 
PAMELA region below ^200 GeV is now d isfavored by the 
new positron measurements with the LAT Ackermann et al.' 
( [201 2b) which indicate that the positron fraction continues to 
rise to this energy. 



It should be noted that the above conclusions are not af- 
fected by the uncertainty in po since both the derived con- 
straints and the region of parameter space compatible with the 
DM interpretation of the CR anomalies scales in the same way 
with pq. The same is true also for the constrains on decaying 
DM. 

Constraints for the case of decaying DM are shown in Fig- 
ure |5] The interpretation of the PAMELA/ferm/ CR features 
in terms of decaying DM is not ruled out in this analysis. The 
limits are stronger than the ones derived in similar analyses 
performed without background modeling (Pa pucci & Strumla] 
2010; Cirelli et al. 2010) and slightly improved over the ones 



derived from observation of Galaxy clusters (Huang et al. 
2012| l. They are compara ble to the limits der i ved from the 



comp arison with the 1GB (Cirelli et al. 2010 Huang et al. 
[20T2| l. 

In Figure [6] we show the CR source distributions derived in 
the fit (the coefficients c"''') and the uncertainty of them. Each 
bin is treated as a parameter for which the profile likelihood 
is built by marginalizing over all the other parameters (linear 
and nonlinear) and DM models and 3 and 5 a uncertainties 
are derived in the usual way. Again, since the Galactic Plane 
and the full energy range is not included in the fit, in inter- 
preting the figure the same caveats as for Figure [3] apply. For 
example, the fitted pCRSD increases in the outerGalaxy: al- 
though protons from the outer Galaxy can propagate to the 
inner Galaxy (while few also propagate from the other side 
of the Galactic Center) and we thus have some sensitivity to 
them in our ROI, clearly, a proper statement on this feature 
would require to fit the outer Galaxy itself. Further caveats 
in interpreting these results include the fact that the CRSDs 
also depend on the ISM distribution and ISRF distribution, 
which are themselves affected by large uncertainties. Finally, 
the CRSDs also depend on the chosen propagation setup. We 
do not attempt here to assess the systematic uncertainties re- 
lated to these issues (see the next section, however, for some 
discussion regarding the ISRF uncertainties). 

Overall, however, it can be noted that the eCRSD and 
pCRSD are in reasonable agreement with each other and with 
typical SNR or pulsars models, except for the rise in the outer 
Galaxy which is not predicted in these models. The plots are 
very similar when no DM contribution is assumed, which is 
expected since DM gives a subdominant contribution to the 
gamma-ray signal 3 kpc away from the GC. 

An interesting point to comment on is how the DM con- 
straints are affected if some of the nuisance parameters, in 
particular the CRSDs, are fixed to benchmark choices taken 
from the literature, instead of being marginalized away. We 
checked the resulting DM constraints for cases in which we 
fix the CRSDs (the eCRSD and /^CRSD are taken to be equal) 



to three different c ommon choices, nam ely SNR (Case & 
Bhattacharya||1998|l, pulsar (jYusifov & Kuguk 2004), aii^ 



simple gaussian model (centered at the Galactic Center, with a 
width of 4.5 kpc) while the rest of the parameters is marginal- 
ized away in the same way as in the main analysis. We find 
that the SNR case gives slightly worse limits (20-30%) while 
pulsar limits are very similar to the results of the main anal- 
ysis (in agreement with the fact that the our best fit CRSDs 
are close to the pulsar distribution), and gaussian limits are 
a factor of <2 better than the pulsar case. This last case 
is understood in light of the fact that a gaussian CRSD, be- 
ing peaked at the Galactic Center, forces the inner Galaxy 
gamma-ray signal to be explained entirely in terms of ordi- 
nary CR sources, leaving little gamma-ray flux to DM and 
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Figure 4. Upper limits on the velocity averaged DM annihilation cross-section including a model of the astrophysical background compared with the limits 
obtained with no modeling of the background. Upper panel: Limits on models in which DM annihilates into bb, for a DM distribution given by the NFW 
distribution (left) and isothermal distribution (right). In the left panel we also add an uncertainty band (red dotted lines) in the 3(t no-background limits which 
would result from varying the local DM density po in the range 0.2-0.7 GeV cm"'. A similar band, not shown in the plot for clarity, would be present for the 
limits including a model of the astrophysical background (see discussion in the text). The horizontal line marks the thermal decoupling cross section expected 
for a generic WIMP candidate. Middle panel: Upper limits for DM annihilation to p^pT . Lower panel: The same, for DM annihilation to t^t^ . The region 
excluded by the analysis with no model of the astrophysical background is indicated in light blue, while the additional region exclu ded by the an alysis with a 
modeling of the background is indicated in light green. The regions of parameter space whic h provide a g ood fit to PAMELA jAdriani et al.|(2009a| (purple) and 
Fermi LAT Abdo et al. (20091 (blue) CR electron and positron data are shown, as derived in |Cirelli et al.|j2010^ and are scaled by a factor of 0.5, to account for 
different assumptions on the local DM density (see text for more details). 
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Figure 5. Lower limits on the lifetime of decaying DM. The panel structure is the same as in Figure[4] In the top left panel we also add an uncertainty band (red 
dotted lines) in the 3(j no-background limits which would result from varying the local DM density po ths range 0.2-0.7 GeV cm"'. A similar band, not shown 
in the plot for clarity, would be present for the limits including a model of the astrophysical background (see discussion in the text). The regions of parameter 
space whi ch provide a good fi t to PAMELA |Adriani et al.|^2009aj l (purple) and Fermi LAT |Abdo et al.|p009j (blue) CR electron and positron data are shown as 
derived in |Ciielli et al.|(2010} and are scaled by a factor oi 1 .4, to account for different assumptions on the local DM density (see text for more details). 
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Figu re 6. Best fit e and p C RSDs and eiTors, marginalized over tlie remaining parameters and over the various DM models considered. The pulsar distribution 
from |Yusifov & Kiigiik|(2004j is also shown for comparison. The source distribution is zeroed within 3 kpc of the Galactic center; see the text. 



thus giving better limits. This is interesting and might indi- 
cate that the DM constraints can become better if indepen- 
dent robust constraints on the CRSD become available. On 
the other hand, as the results of the main analysis show, the 
fitted CRSD does not favor a Gaussian distribution (at least 
at a qualitative level; we have not performed a quantitative 
comparison for the reasons explained above). For the time 
being, thus, the use of CRSD as a nuisance parameter seems 
the best approach, which, at the same time, leaves freedom in 
the fit to explore the degeneracy with DM and limit the fit to 
explore CRSDs which are in reasonable agreement with the 
(gamma-ray) data. 

Finally, we show in Figure [7] the counts map in our ROI, 
together with the model prediction for a model close to the 
best fit {zh = 10, kpc 7^,2 = 2.3 and d2Hl=0.0140 x lO'^" mag 
cm^) and its residuals when DM (of mass m^=15Q GeV and 
annihilating into bb) is included or not in the fit. It can be seen 
that the residuals are mostly flat, meaning that the model (and 
the models close to the minimum) are a reasonable fit to the 
data. A few features are however present, like the excess in 
the vicinity of (Z,/7)~(-45,10) and ~(7,-15) which seem to be 
related to the low latitude tip of Loop 1 and to the low latitude 
part of the South Lobe/Bubble, respectively. The two promi- 
nent negative residuals near (I,b)^(-15,5) and ~(20,-10), in- 
stead, approximately contour the Lobes and thus seem to be 
an artifact of the fit to compensate for this missing component. 
Gas misplaced in incorrect annuli also could be an alternative 
explanation. 

We also show the point-source mask used based on the 
IFGL catalog and, for co mparison, the mask based on the 
2FGL catalog ( Nolan et al. 20 12) and the residuals using this 
mask. Overall, it can be seen that the 2FGL mask covers 
few point sources which are apparent in the residuals with 
the IFGL mask. The large scale features in the residuals are 
however unchanged, apart from a small part of Loop 1 near 
(Z,/7)~(-45,10) which is resolved into sources. 

9. DISCUSSION ON MODEL UNCERTAINTIES 

In deriving our limits above we have taken into account 
many possible uncertainties like the ones in the e,/?CRSDs, 
in Zh, the electron index and the dust to gas ratio. We check 
below the importance of further uncertainties which we have 
not considered expHcitly in our scan. 



An important component for which there is still a consid- 
erable uncertainty is the ISRF. In particular, the ISRF in the 
inner Galaxy is quite uncertain and the default model we used 
could be a substantial underestimate of the true one in this 
region. Very different ISRFs would affect the propagation of 
CREs through energy losses and this could be especially rele- 
vant for the DM models in which the IC component is impor- 
tant and provide strong constraints, like /i channels. Modes 
dominated by prompt radiation, like b and t should, instead, 
not be significantly subject to uncertainties in the ISRF. To 
make an explicit check we repeated ou r entire analysis us- 
ing a different ISRF model ( |Porter|2012| l, which has the bulge 
component increased by a factor of 10 (see [Forter et al.| ( |2008| l 
for a detailed definition), which implies an overall increase 
in the inner Galaxy of a factor of 2. The DM limits with 
this enhanced ISRF were, however, not appreciably affected. 
We verified that the enhanced ISRF produces an enhanced 
IC component, but only within a few degrees of the Galactic 
Center, thus not affecting the fit in our ROI. It also should be 
stressed that a more intense ISRF implies more IC emission 
for the DM IC too, so that assuming a lower ISRF gives con- 
servative limits. Finally, an ISRF lower than the one assumed 
here is also p ossible, as the results obtained in |Ackermann] 
[eraLld20T2a| l (see Figure 11 there) for the CRSD following 
the pulsar distribution seems to indicate. However, the "ISRF 
normalization" reported in Ackermann et al. (2012a), is more 
precisely a proxy for a combination of ISRF intensity, normal- 
ization of the CRE spectrum and halo size, so that alternative 
explanations are possible. 

We also checked more systematically other sources of un- 
certainties, but in a more simplified setup: we set a partic- 
ular model as reference and then we varied each parameter 
one at a time, keeping the others fixed, and for each case we 
calculate the percentage variation in DM limits for selected 
DM models. We vary the parameters derived from the CR fit, 
va, Jpj, lp.2, Phr.p, and t he (DQ,Zh) relation , within the uncer- 
tainty ranges derived in l Ackermann et al.| ( |2012a| l enlarging 
it by a factor of ~2 to take mto account possible systematic 



uncertainties (the errors quoted in Ackermann et al. ( 2012a| 
are statistical only). We also include in the list of the tested 
parameters the ones which are included in our model scan 
(CRSD, d2HI, 7f 2, {Do,Zh)) to allow for a direct comparison. 
The following set of parameters, which lie close to the best fit 
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Figure 7. Counts map of the ROI that we consider (upper left panel), model prediction for a model (without DM) close to the best fit (zi, = 10 kpc, 7^ 2 = 2.3 
and d2HI=0.0140 X 10"^" mag cm") parameter region (upper right), and residuals in units of a for the same model (second row left) and when DM (of mass 
m-^=l50 GeV and annihilating into hb) is also included in the fit (second row right). Third row: same as second row but with 2FGL point sources masked instead 
of IFGL. Fourth row: IFGL mask (left) and 2FGL mask (right). The model and data counts and the residuals have been smoothed with a 1.25° Gaussian filter. 
The point sources mask in the residuals have been applied before and after the smoothing. 



of our analysis, was chosen for the reference model: =36 
km s-i. Do = 5.0 lO^^cm^s-i, = 4 kpc, 6 = 0.3, 7^,1 = 1.9, 
7p,2 = 2.39, ptr.p = 1 1.5, 7,,2 = 2.45, d2HI=0.014 x lO'^" mag 
cm^, CRSD=SNR, K=0 km s"^ Results ai-e shown in Table 

m 

We can see that CR parameters such as va and 7^.1, 7^,2, 
Pbr.p, Vc and even different gas maps have very low (< 10%) 
impact on the DM limits. The table confirms that -fe,! and 
the CRSD (which we fix here to be the same for protons and 
CREs) are the main parameters degenerate with DM and thus 
affecting the limits the most (up to 60%). The diffusion con- 
stant Do is tightly correlated to the halo height z/,. Therefore 
we vary the parameter pair (Do, z/,) instead of the single pa- 
rameters, using their relation derived from the fit to the CR 
data described in Sec. |2] Nonetheless, the combination Do 
and z/, are included individually in the parameter scans of 
the previous section used for the main results. As an addi- 
tional check of the effect of the CR propagation parameters 
on the DM limits, we find DM limits in three theoretical CR 
propagation setups: plain diffusion (PD, characterized by in- 
dex of diffusion of 5 = 0.6), Kraichanian (KRA, S = 0.5), and 
Kolmogorov (KOL, 6 = 0.3). In these cases, the rest of the 
CR propagation parameters are found from the best fit to the 



(|2012a|;[d rBernardo etaLl pOTT]!. These fits to CR are per- 
formed again without any DM component. We find that DM 
limits in these three CR diffusion setups are also barely af- 
fected, in particular when compared to the effect of the CR 
source distribution, as shown in Table [3] 

Finally we also consider an alternative configuration of the 
Galactic Magnetic Field (GMF). The reference one is the de- 
fault configuration used in GALPROP with an exponential 
profile in R and z and length scales of 10 kpc in R and 2 
kpc in z, normalized to 5 /iG locally (Conf 1). The alterna- 
tive configuration we tested has in addition a further compo- 
nent of constant 100 ^G intensity within 0.4 kpc from the 



CR data following the method described in Ackermann et al. 



Galactic Center, as motivated by a recent work Crocker et al. 
(2010) (Conf 2). This alternative configuration also produces 
changes in the limits of less than 10%. 

10. SUMMARY 

In this work we constrain the contribution to the diffuse 
gamma-ray emission from DM annihilating or decaying in the 
Milky Way halo, based on the Fermi LAT gamma-ray data. 
We first present the most conservative limits on DM assum- 
ing that all LAT photons from the Halo are produced by an- 
nihilating/decaying WlMPs. Then, based on our current best 
knowledge of the Galactic diffuse emission (Ackermann et aLj 
|2012a| l, we use GALPROP to model the astrophysical diffuse 
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Parameter 


\5ct/ct\ [%lhh 


\5a/a\ [%], 


I'A [ 30; 36; 45] km s"' 

7^,1 [1.8; 1.9; 2;] 
7^,2 [2.35; 2.39; 2.45] 
Pbr.p [ 10; 11.5; 12.5] GV 
d2HI [ 0.01 10, 0.0140; 0.0170] IQ-^o mag cm- 


[6; 0; 11] 
[ 1.0; 0; 2.5] 
[2.5;0; 1.5] 
[0.5;0; 1.0] 

[3;0; 12] 


[4.;0;9] 
[1.5; 0; 2.0] 
[2.5; 0; 1.5] 
[0.9; 0; 1.5] 

[ 3;0; 9] 


7^.2 [ 2.0; 2.45; 2.6] 
(Do.zh) [ (5.0e28, 4); (7.1e28, 10)] cm-s"' 
CRSD [ SNR; Pulsar] 


[ 17; 0; 7] 
[ 0; 10] 
[0;61] 


[ 18; 0; 5] 
[0; 7] 
[ 0; 59] 


KRA{S = 0.5); KOL((5 = 0.3); PD(5 = 0.6) 
Vc [0; 20] km s"' 
GMF [ Conf 1, Conf 2] 


[ 4.0; 0; 3.0] 
[0; 6] 
[0; 3] 


[1.0; 0; 5] 
[0; 4] 
[0; 8] 



Table 3 

Relative variation |(5(t/(t|[%] of the limits on the DM velocity averaged annihilation cross-section derived in this work with respect to changes in the underlying 
astrophysical diffuse emission model. The table shows the relative variation for selected DM models (bh and fi'^fi" channel, for a 150 GeV DM) in a simplified 
set-up when only one parameter is varied at a time. Each row corresponds to the indicated parameter. The bold values correspond to the reference value. 



background, and, using a profile likelihood approach, we ex- 
plore the effects of various poorly constrained parameters in 
the modeling of the astrophysical background, e.g., the dif- 
fusive halo height, the CR source distribution, the index of 
the electron injection spectrum and the dust to gas ratio in 
order to get more robust constraints. We also remove astro- 
physical CR sources within 3 kpc from the Galactic Center so 
that any potential astrophysical contribution in this region is 
attributed to DM, resulting in more conservative constraints. 
Overall, rather than being due to residual astrophysical model 
uncertainties, the remaining major uncertainties in the DM 
constraints from the Halo region come from the modeling of 
the DM signal itself. The main uncertainty is in the normaliza- 
tion of the DM profile, which is fixed through the local value 
of the DM densi ty. We use the recent determination po = 0.43 
GeV cm"^ from jSalucci et al. ( 2010 1, which has, however, a 
large uncertainty, with values in the range 0.2-0.7 GeV cm"^ 
still viable. A large uncertainty in po is particular important 
for annihilation constraints since they scale like pg, while for 
constraints on decaying DM the scaling is only linear A less 
important role is played by the uncertainties in the DM profile, 
since in the Halo region different profiles predict similar DM 
densities. When using the lowest allowed DM density used in 
literature of po = 0.2 GeV cm"-' (see e.g. Salucci et al. 2010 1, 
our limits worsen by a factor 4 (2) for annihilating (decaying) 
DM, as shown in Figures H] and l5] A better determination of 
the local DM density, as well as or the parameters determining 
the global structure of the DM Halo, is therefore of the utmost 
importance for reducing the uncertainties related to DM con- 
straints from DM halo, but it is beyond the scope of this paper 
and is the subject of dedicated studies. 

Bearing this in mind, the limits we obtain are competitive 
with complementary probes of DM like dwarfs, clusters of 
galaxies, or recent constraints obtained from CMB observa- 
tions for DM models with prompt spectra, and significantly 
improve over these studies for DM models with significant IC 
contribution such as DM annihilating into The limits 

we derive for leptonic models challenge the interpretation of 
the RA^MELA and Fermi CR anomalies as annihilation of DM 
in the Galactic Halo, while they are not constraining enough 
to exclude the interpretation in terms of decaying DM. We 
note that this last conclusion is not affected by the uncertainty 
in Po since both the derived constraints and the region of pa- 
rameter space compatible with the DM interpretation of the 
CR anomalies scale in the same way with po. 



An obvious improvement of this analysis would be a full 
scan of the CR parameters from a simultaneous fit to both 
gamma and cosmic-ray data. An effort in this direction is 
currently ongoing and will be reported in future publications. 
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